Piecewise Linear Functions Research Articles

Data-Driven Initial Gap Identification of Piecewise-Linear Systems Using Sparse Regression and Universal Approximation Theorem

Abstract This paper proposes a method for identifying an initial gap in piecewise-linear systems from data. Piecewise-linear systems appear in many engineered systems such as degraded mechanical systems and infrastructures, and are known to show strong nonlinearities. To analyze the behavior of such piecewise-linear systems, it is necessary to identify the initial gap, at which the system behavior switches. The proposed method identifies the initial gap by discovering the governing equations using sparse regression and calculating the gap based on the universal approximation theorem. A key step to achieve this is to approximate a piecewise-linear function by a finite sum of piecewise-linear functions in sparse regression. Equivalent gap is then calculated from the coefficients of the multiple piecewise-linear functions and their respective switching points in the obtained equation. The proposed method was first applied to a numerical model to confirm its applicability to piecewise-linear systems. Experimental validation of the proposed method has then been conducted with a simple mass-spring-hopping system, where the method successfully identified the initial gap in the system with high accuracy.

Bifurcation structures of a two-dimensional piecewise linear discontinuous map: analysis of a cobweb model with regime-switching expectations

AbstractWe consider the bifurcations occurring in a two-dimensional piecewise-linear discontinuous map that describes the dynamics of a cobweb model in which firms rely on a regime-switching expectation rule. In three different partitions of the phase plane, separated by two discontinuity lines, the map is defined by linear functions with the same Jacobian matrix, having two real eigenvalues, one of which is negative and one equal to 0. This leads to asymptotic dynamics that can belong to two or three critical lines. We show that when the basic fixed point is attracting, it may coexist with at most three attracting cycles. We have determined their existence regions, in the two-dimensional parameter plane, bounded by border collision bifurcation curves. At parameter values for which the basic fixed point is repelling, chaotic attractors may exist - either one that is symmetric with respect to the basic fixed point, or, if not symmetric, the symmetric one also exists. The homoclinic bifurcations of repelling cycles leading to the merging of chaotic attractors are commented by using the first return map on a suitable line. Moreover, four different kinds of homoclinic bifurcations of a saddle 2-cycle, leading to divergence of the generic trajectory, are determined.

A new flexible and partially monotonic discrete choice model

The poor predictability and the misspecification arising from hand-crafted utility functions are common issues in theory-driven discrete choice models (DCMs). Data-driven DCMs improve predictability through flexible utility specifications, but they do not address the misspecification issue and provide untrustworthy behavioral interpretations (e.g., biased willingness to pay estimates). Improving interpretability at the minimum loss of flexibility/predictability is the main challenge in the data-driven DCM. To this end, this study proposes a flexible and partially monotonic DCM by specifying the systematic utility using the Lattice networks (i.e., DCM-LN). DCM-LN ensures the monotonicity of the utility function relative to the selected attributes while learning attribute-specific non-linear effects through piecewise linear functions and interaction effects using multilinear interpolations in a data-driven manner. Partial monotonicity could be viewed as domain-knowledge-based regularization to prevent overfitting, consequently avoiding incorrect signs of the attribute effects. The light architecture and an automated process to write monotonicity constraints make DCM-LN scalable and translatable to practice. The proposed DCM-LN is benchmarked against deep neural network-based DCM (i.e., DCM-DNN) and a DCM with a hand-crafted utility in a simulation study. While DCM-DNN marginally outperforms DCM-LN in predictability, DCM-LN highly outperforms all considered models in interpretability, i.e., recovering willingness to pay at individual and population levels. The empirical study verifies the balanced interpretability and predictability of DCM-LN. With superior interpretability and high predictability, DCM-LN lays out new pathways to harmonize the theory-driven and data-driven paradigms.

The bifurcation structure within robust chaos for two-dimensional piecewise-linear maps

We study two-dimensional, two-piece, piecewise-linear maps having two saddle fixed points. Such maps reduce to a four-parameter family and are well known to have a chaotic attractor throughout open regions of parameter space. The purpose of this paper is to determine where and how this attractor undergoes bifurcations. We explore the bifurcation structure numerically by using Eckstein’s greatest common divisor algorithm to estimate from sample orbits the number of connected components in the attractor. Where the map is orientation-preserving the numerical results agree with formal results obtained previously through renormalisation. Where the map is orientation-reversing or non-invertible the same renormalisation scheme appears to generate the bifurcation boundaries, but here we need to account for the possibility of some stable low-period solutions. Also the attractor can be destroyed in novel heteroclinic bifurcations (boundary crises) that do not correspond to simple algebraic constraints on the parameters. Overall the results reveal a broadly similar component-doubling bifurcation structure in the orientation-reversing and non-invertible settings, but with some additional complexities.

Solution path algorithm for distributionally robust regression

ABSTRACT In this paper, we propose a general distributionally robust regression model based on distributionally robust optimization theory. The proposed model has a piecewise linear loss function and elastic net penalty term, and it generalizes many other regression models. We prove the piecewise linear property of the optimal solutions to this model, which enables us to develop a solution path algorithm for the hyperparameter tuning. A Doubly regularized Least Absolute Deviations (DrLAD) regression model is proposed based on this framework, and a solution path algorithm is developed to speed up the tuning of two hyperparameters in this model. Numerical experiments are implemented to validate the performance of this model and the computational efficiency of the solution path algorithm.

Gauge-Optimal Approximate Learning for Small Data Classification.

Small data learning problems are characterized by a significant discrepancy between the limited number of response variable observations and the large feature space dimension. In this setting, the common learning tools struggle to identify the features important for the classification task from those that bear no relevant information and cannot derive an appropriate learning rule that allows discriminating among different classes. As a potential solution to this problem, here we exploit the idea of reducing and rotating the feature space in a lower-dimensional gauge and propose the gauge-optimal approximate learning (GOAL) algorithm, which provides an analytically tractable joint solution to the dimension reduction, feature segmentation, and classification problems for small data learning problems. We prove that the optimal solution of the GOAL algorithm consists in piecewise-linear functions in the Euclidean space and that it can be approximated through a monotonically convergent algorithm that presents-under the assumption of a discrete segmentation of the feature space-a closed-form solution for each optimization substep and an overall linear iteration cost scaling. The GOAL algorithm has been compared to other state-of-the-art machine learning tools on both synthetic data and challenging real-world applications from climate science and bioinformatics (i.e., prediction of the El Niño Southern Oscillation and inference of epigenetically induced gene-activity networks from limited experimental data). The experimental results show that the proposed algorithm outperforms the reported best competitors for these problems in both learning performance and computational cost.

Explicit machine learning-based model predictive control of nonlinear processes via multi-parametric programming

Machine learning-based model predictive control (ML-MPC) has been developed to control nonlinear processes with unknown first-principles models. While ML models can capture nonlinear dynamics of complex systems, the complexity of ML models leads to increased computation time for real-time implementation of ML-MPC. To address this issue, in this work, we propose an explicit ML-MPC framework for nonlinear processes using multi-parametric programming. Specifically, a self-adaptive approximation algorithm is first developed to obtain a piecewise linear affine function that approximates the behaviors of ML models. Then, multi-parametric quadratic programming (mpQP) problems are formulated to generate the solution map for states in discretized state–space. Furthermore, to accelerate the implementation of explicit ML-MPC, a neighbor-first search algorithm is developed. Finally, an example of a chemical reactor is used to demonstrate the effectiveness of the explicit ML-MPC.

Eye movement analysis for real-world settings using segmented linear regression

Eye movement analysis is critical to studying human brain phenomena such as perception, cognition, and behavior. However, under uncontrolled real-world settings, the recorded gaze coordinates (commonly used to track eye movements) are typically noisy and make it difficult to track change in the state of each phenomenon precisely, primarily because the expected change is usually a slower transient process. This paper proposes an approach, Improved Naive Segmented linear regression (INSLR), which approximates the gaze coordinates with a piecewise linear function (PLF) referred to as a hypothesis. INSLR improves the existing NSLR approach by employing a hypotheses clustering algorithm, which redefines the final hypothesis estimation in two steps: (1) At each time-stamp, measure the likelihood of each hypothesis in the candidate list of hypotheses by using the least square fit score and its distance from the k-means of the hypotheses in the list. (2) Filter hypothesis based on a pre-defined threshold. We demonstrate the significance of the INSLR method in addressing the challenges of uncontrolled real-world settings such as gaze denoising and minimizing gaze prediction errors from cost-effective devices like webcams. Experiment results show INSLR consistently outperforms the baseline NSLR in denoising noisy signals from three eye movement datasets and minimizes the error in gaze prediction from a low precision device for 71.1% samples. Furthermore, this improvement in denoising quality is further validated by the improved accuracy of the oculomotor event classifier called NSLR-HMM and enhanced sensitivity in detecting variations in attention induced by distractor during online lecture.

The hidden sensitivity of non-smooth dynamics

Switches in dynamical systems are known to exhibit wildly different behaviours depending on how they are modelled, for instance whether they occur as smooth transitions or abrupt jumps, and whether they involve delays or discrete perturbations. These differences arise because switches are sensitive to perturbation, but there is limited knowledge about where this sensitivity comes from. Here we take a switch in a simple one-dimensional system, then discretise time and introduce a small delay. The resulting system is described by a piecewise-linear map with incredibly complex dynamics, and is sensitive to parameter changes even if the time-steps and delays are made infinitesimally small. We show that this sensitivity reveals itself in the more versatile numerical tool of the transition count, which captures the likelihood of switching occurring at any instant in a simulation. We use this to show that sensitivity to parameters persists in a system with two switches, where it then has large-scale dynamical effects. The use of transition counts in this way may prove a versatile tool for studying more complex switching processes in general.

Robust space-time finite element methods for parabolic distributed optimal control problems with energy regularization

As in our previous work (SINUM 59(2):660–674, 2021) we consider space-time tracking optimal control problems for linear parabolic initial boundary value problems that are given in the space-time cylinder Q=Ω×(0,T)\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$Q = \\Omega \ imes (0,T)$$\\end{document}, and that are controlled by the right-hand side zϱ\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$z_\\varrho $$\\end{document} from the Bochner space L2(0,T;H-1(Ω))\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$L^2(0,T;H^{-1}(\\Omega ))$$\\end{document}. So it is natural to replace the usual L2(Q)\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$L^2(Q)$$\\end{document} norm regularization by the energy regularization in the L2(0,T;H-1(Ω))\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$L^2(0,T;H^{-1}(\\Omega ))$$\\end{document} norm. We derive new a priori estimates for the error ‖u~ϱh-u¯‖L2(Q)\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\\Vert \\widetilde{u}_{\\varrho h} - \\overline{u}\\Vert _{L^2(Q)}$$\\end{document} between the computed state u~ϱh\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\\widetilde{u}_{\\varrho h}$$\\end{document} and the desired state u¯\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\\overline{u}$$\\end{document} in terms of the regularization parameter ϱ\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\\varrho $$\\end{document} and the space-time finite element mesh size h, and depending on the regularity of the desired state u¯\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\\overline{u}$$\\end{document}. These new estimates lead to the optimal choice ϱ=h2\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\\varrho = h^2$$\\end{document}. The approximate state u~ϱh\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\\widetilde{u}_{\\varrho h}$$\\end{document} is computed by means of a space-time finite element method using piecewise linear and continuous basis functions on completely unstructured simplicial meshes for Q. The theoretical results are quantitatively illustrated by a series of numerical examples in two and three space dimensions. We also provide performance studies for different solvers.

Spurious Local Minima are Common for Deep Neural Networks With Piecewise Linear Activations.

In this article, theoretically, it is shown that spurious local minima are common for deep fully connected networks and average-pooling convolutional neural networks (CNNs) with piecewise linear activations and datasets that cannot be fit by linear models. Motivating examples are given to explain why spurious local minima exist: each output neuron of deep fully connected networks and CNNs with piecewise linear activations produces a continuous piecewise linear (CPWL) function, and different pieces of the CPWL output can optimally fit disjoint groups of data samples when minimizing the empirical risk. Fitting data samples with different CPWL functions usually results in different levels of empirical risk, leading to the prevalence of spurious local minima. The results are proved in general settings with arbitrary continuous loss functions and general piecewise linear activations. The main proof technique is to represent a CPWL function as maximization over minimization of linear pieces. Deep networks with piecewise linear activations are then constructed to produce these linear pieces and implement the maximization over minimization operation.

Fixed-Time Adaptive Parameter Estimation for Hammerstein Systems Subject to Dead-Zone

The Hammerstein model has proven its ability for modeling many industrial servo systems, but the corresponding parameter estimation remains as a challenging problem, especially for the purpose of achieving fast convergence. This paper presents a fixed-time (FxT) adaptive parameter estimation scheme for Hammerstein systems with an asymmetric dead-zone to ensure the convergence of estimated parameters in fixed-time. A continuous piecewise linear (CPL) function is adopted to model the nonlinear dead-zone dynamics to remedy the use of intermediate variable. To avoid using the system state information, a K-filter is used to construct a relationship between the unknown system parameters and the input/output signals. Then a new adaptive law based on the parameter error and the FxT scheme is developed to online estimate the lumped unknown parameters and ensure the fixed-time convergence, where an online updated auxiliary variable of regressor is suggested to eliminate the impact of the regressor on the convergence performance. Both numerical simulations and experimental results are given to demonstrate the effectiveness of proposed methods.

Analysis of dynamics of a map-based neuron model via Lorenz maps.

Modeling nerve cells can facilitate formulating hypotheses about their real behavior and improve understanding of their functioning. In this paper, we study a discrete neuron model introduced by Courbage et al. [Chaos 17, 043109 (2007)], where the originally piecewise linear function defining voltage dynamics is replaced by a cubic polynomial, with an additional parameter responsible for varying the slope. Showing that on a large subset of the multidimensional parameter space, the return map of the voltage dynamics is an expanding Lorenz map, we analyze both chaotic and periodic behavior of the system and describe the complexity of spiking patterns fired by a neuron. This is achieved by using and extending some results from the theory of Lorenz-like and expanding Lorenz mappings.

Optimal Mechanism in a Dynamic Stochastic Knapsack Environment

This study introduces an optimal mechanism in a dynamic stochastic knapsack environment. The model features a single seller who has a fixed quantity of a perfectly divisible item. Impatient buyers with a piece-wise linear utility function arrive randomly and they report the two-dimensional private information: marginal value and demanded quantity. We derive a revenue-maximizing dynamic mechanism in a finite discrete time framework that satisfies incentive compatibility, individual rationality, and feasibility conditions. This is achieved by characterizing buyers' utility and utilizing the Bellman equation. Moreover, we establish the essential penalty scheme for incentive compatibility, as well as the allocation and payment policies. Lastly, we propose algorithms to approximate the optimal policy, based on the Monte Carlo simulation-based regression method and reinforcement learning.

Piecewise Linear Transformation – Propagating Aleatoric Uncertainty in Neural Networks

Real-world data typically exhibit aleatoric uncertainty which has to be considered during data-driven decision-making to assess the confidence of the decision provided by machine learning models. To propagate aleatoric uncertainty represented by probability distributions (PDs) through neural networks (NNs), both sampling-based and function approximation-based methods have been proposed. However, these methods suffer from significant approximation errors and are not able to accurately represent predictive uncertainty in the NN output. In this paper, we present a novel method, Piecewise Linear Transformation (PLT), for propagating PDs through NNs with piecewise linear activation functions (e.g., ReLU NNs). PLT does not require sampling or specific assumptions about the PDs. Instead, it harnesses the piecewise linear structure of such NNs to determine the propagated PD in the output space. In this way, PLT supports the accurate quantification of predictive uncertainty based on the criterion exactness of the propagated PD. We assess this exactness in theory by showing error bounds for our propagated PD. Further, our experimental evaluation validates that PLT outperforms competing methods on publicly available real-world classification and regression datasets regarding exactness. Thus, the PDs propagated by PLT allow to assess the uncertainty of the provided decisions, offering valuable support.

Vibration Suppression of Piecewise-Linear Stiffness Nes System Under Random Excitation

Nonlinear energy sinks (NESs) are critically important for structural vibration suppression. They can absorb vibrational energy across a broad frequency spectrum, possess strong robustness, and have a relatively small mass. This study addresses the vibration suppression in a piecewise linear stiffness NES system under random excitation. Initially, a theoretical model of the piecewise linear stiffness NES system is developed. The piecewise linear stiffness function is approximated using Legendre polynomial approximation. Following this, the steady-state Fokker–Planck–Kolmogorov (FPK) equation of the system is formulated via the Generalized Harmonic Function Method. The FPK equation is solved using the fourth-order central difference method, and the effectiveness of this approach is validated by comparing the FDM results with numerical simulations. Lastly, the influence of varying system parameters on the stability of the piecewise linear stiffness NES system is analyzed.

Analytical studies on third-order chaotic systems with Sprott type nonlinearities and their microcontroller implementation

The evolution of chaos in a generic third-order autonomous mathematical model with nonlinearities described by simple mathematical functions is reported in this paper. The nonlinearities termed as Sprott type nonlinear functions are used in the design of a class of third-order systems exhibiting chaotic behavior. The evolution and confirmation of chaos in their system dynamics is observed through numerical simulation studies of one-parameter bifurcation diagrams and Lyapunov exponents. Analytical solutions are developed for systems with piecewise-linear nonlinear functions. Finally, the microcontroller implementation of the third-order system equations with different nonlinearities and analog circuit simulation results are presented to confirm the numerical and analytical results. Chaos in generic third-order systems studied through numerical, analytical and microcontroller results has been reported in the literature for the first time.

The necessity of the sausage-string structure for mode-locking regions of piecewise-linear maps

Piecewise-smooth maps are used as discrete-time models of dynamical systems whose evolution is governed by different equations under different conditions (e.g. switched control systems). By assigning a symbol to each region of phase space where the map is smooth, any period-p solution of the map can be associated to an itinerary of p symbols. As parameters of the map are varied, changes to this itinerary occur at border-collision bifurcations (BCBs) where one point of the periodic solution collides with a region boundary. It is well known that BCBs conform broadly to two cases: persistence, where the symbolic itinerary of a periodic solution changes by one symbol, and a nonsmooth-fold, where two solutions differing by one symbol collide and annihilate. This paper derives new properties of periodic solutions of piecewise-linear continuous maps on Rn to show that under mild conditions BCBs of mode-locked solutions on invariant circles must be nonsmooth-folds. This explains why Arnold tongues of piecewise-linear maps exhibit a sausage-string structure whereby changes to symbolic itineraries occur at codimension-two pinch points instead of codimension-one persistence-type BCBs. But the main result is based on the combinatorical properties of the itineraries, so the impossibility of persistence-type BCBs also holds when the periodic solution is unstable or there is no invariant circle.

Optimal scheduling of an integrated energy system considering carbon trading mechanism and multi-energy supply uncertainties in a real-time price environment

This research presents an optimized scheduling model for Integrated Energy Systems (IES) that harmonizes economic efficiency with low carbon emissions in a Real-Time Price (RTP) environment. Integrating a ladder-type carbon tax approach, this model adeptly addresses the unpredictability of renewable energy generation and IES load demands. It introduces a carbon trading mechanism for efficient emission control, complemented by an RTP model that manages interactions between the IES and the main power grid through piecewise linear functions. The model also proposes an adaptable operation strategy for thermoelectric and electric cooling ratios, enhancing the IES’s economic and operational flexibility. The objective function, crafted on a deviation preference optimization basis, reflects the diverse goals of decision-makers, and the model is skillfully translated into a mixed-integer linear programming format through a series of mathematical transformations. Simulation results demonstrate the model’s efficacy in achieving an optimal balance between cost-effectiveness and reliability, particularly by setting a strategic confidence level for spinning reserves. Additionally, incorporating a ladder-type carbon tax and Electric Vehicle (EV) charging within the IES underscores its significant environmental benefits, marking a stride towards a low-carbon future.

Chance-Constrained Optimization for a Green Multimodal Routing Problem with Soft Time Window under Twofold Uncertainty

This study investigates a green multimodal routing problem with soft time window. The objective of routing is to minimize the total costs of accomplishing the multimodal transportation of a batch of goods. To improve the feasibility of optimization, this study formulates the routing problem in an uncertain environment where the capacities and carbon emission factors of the travel process and the transfer process in the multimodal network are considered fuzzy. Taking triangular fuzzy numbers to describe the uncertainty, this study proposes a fuzzy nonlinear programming model to deal with the specific routing problem. To make the problem solvable, this study adopts the fuzzy chance-constrained programming approach based on the possibility measure to remove the fuzziness of the proposed model. Furthermore, we use linear inequality constraints to reformulate the nonlinear equality constraints represented by the continuous piecewise linear functions and realize the linearization of the nonlinear programming model to improve the computational efficiency of problem solving. After model processing, we can utilize mathematical programming software to run exact solution algorithms to solve the specific routing problem. A numerical experiment is given to show the feasibility of the proposed model. The sensitivity analysis of the numerical experiment further clarifies how improving the confidence level of the chance constraints to enhance the possibility that the multimodal route planned in advance satisfies the real-time capacity constraint in the actual transportation, i.e., the reliability of the routing, increases both the total costs and carbon emissions of the route. The numerical experiment also finds that charging carbon emissions is not absolutely effective in emission reduction. In this condition, bi-objective analysis indicates the conflicting relationship between lowering transportation activity costs and reducing carbon emissions in routing optimization. The sensitivity of the Pareto solutions concerning the confidence level reveals that reliability, economy, and environmental sustainability are in conflict with each other. Based on the findings of this study, the customer and the multimodal transport operator can organize efficient multimodal transportation, balancing the above objectives using the proposed model.